Introduction
Soil salinization is not only a leading cause of the
deteriorating ecological environment, but also a major abiotic factor affecting
crop yield around the world (Zhu 2001). According to related statistics, there
is 4.0×108 hm2 saline-alkali land in the world, of which,
an area of 3.6×107 hm2 is located in China (Zhang et al. 2007). The cultivation of new
salt-tolerance plant varieties is an effective way to utilize the saline-alkali
land. Meanwhile, investigating the salt-tolerance mechanisms in plant and
identifying the salt-tolerance genes can lay a solid foundation for cultivating
the new salt-tolerance plant varieties.
Research on the salt-tolerance
mechanisms of plants had been carried out over the past many decades. Selective
ion absorption and compartmentalization play key roles in maintaining ion
homeostasis in cytoplasm. For instances, salt overly sensitive (SOS) and Na+/H+
antiporter (NHX) can keep lower Na+ content in cytoplasm (Zhu 2003),
whereas the K+ transporters with high affinity (HKTs) can improve K+
and limit Na+ transportation from root to leaf (Tang et al. 2015). Salt stress can induce ROS
production, including hydrogen peroxide (H2O2),
superoxide radicals (•O2−) and hydroxyl radicals
(•OH), eventually causing oxidative damage to cytomembrane (Miller et al. 2008). For alleviating the
ROS-induced peroxidation damage, two kinds of ROS scavenging systems have
evolved in plant, including the non-enzymatic antioxidants and the enzymatic
antioxidants. The non-enzymatic antioxidants include glutathione (GSH) and
ascorbic acid (AsA) etc.; whereas the enzymatic
antioxidants consist of peroxidase (POD), catalase (CAT) and superoxide
dismutase (SOD) etc., and they play crucial roles in altering the ROS
homeostasis (Deinlein et al. 2014).
Salt stress leads to imbalanced osmotic regulation; as a result, some
osmolytes, such as free proline, sugar and betaine, are synthesized in cells to
regulate the osmotic balance in plants (Ingram and Bartels 1996; Ashraf and
Foolad 2007). In addition, the expression of salt-tolerance genes in plants are
regulated by transcription factors (TFs), and many of which, including DREB,
MYB, AP2/ERF and NAC families etc, exert vital parts in plant tolerance to salt
stress (Deinlein et al. 2014).
Halophyte is a plant that grows
regularly and completes the life cycles under the single salt concentration of
>70 mmol·L-1 (Flowers and Colmer 2008). Z. macrostachya, the perennial warm-season turfgrass native to
China, Japan as well as Korean Peninsula, and mainly grows in the coastal wetlands
of Shandong, Jiangsu and Zhejiang provinces in China. Z. macrostachya can rapidly spread through rhizomes and stolons to
form the dense turf with the deep root system, which can thereby be used as the
soil-conserving, dike-protecting and sand-fixing turf. Moreover, Z. macrostachya is an euhalophyte, which
exhibits tolerance to salinity and may be potentially used in landscaping of
saline-alkali land. According to our previous research, Z. macrostachya tolerated 355 mM NaCl
stress (Hu and Zhang 2010) and further research found
that the roots of Z. macrostachya
limited Na+ absorption while improved K+ absorption and
transportation from roots to leaves, accumulated free proline and soluble
sugars and enhanced the POD activity under salt stress, and all of these
improved the salt-tolerance of Z.
macrostachya (Hu and Zhang 2009, 2010; Hu et al. 2016). Nonetheless, the above studies are limited to
physiological indexes, and the molecular mechanisms of salt-tolerance in this
species remain unclear so far.
RNA-seq has emerged as a powerful
tool to analyze genes expressional changes, which reflects the molecular
mechanisms of plant response to salt stress, and has been utilized for investigating
the molecular mechanisms of salt-tolerance in many halophytes, such as Iris lactea var. Chinensis (Gu et al. 2018), Prunellae Spica (Liu et al.
2020) and Rhizophora mucronata (Meera
and Augustine 2020). To date, changes in global genes expression of Z. macrostachya under salt stress are
still unknown, which limits the understanding towards the molecular mechanisms
of salt-tolerance in this species. In this study, the Illumina HiSeq XTen
sequencing platform was used to generate a roots reference transcriptome
dataset and to explore DEGs with aims to improve our comprehensive
understanding of the mechanisms of salt-tolerance at transcriptome level and to
identify the DEGs involved in the salt-tolerance of Z. macrostachya.
Plant material, salt stress and RNA extraction
Zoysia
macrostachya samples were collected from coastal wetland located 45.0
km east of Yancheng of Jiangsu province, China. Samples were brought back to
the laboratory and planted in six PVC tubes (35 cm in length and 15 cm in
diameter) filled with river sand, respectively. The plants were grown in a
growth chamber with a 12 h light/12 h dark cycle, 30/20℃ day/night
temperature, 800 μmol m-2·s-1 light
intensity and a relative humidity of 80%. The plants were watered at intervals
of three days, and irrigated weekly with 200 mL 1/2 Hoagland's nutrient
solution during growth. After 30 d of growth, the plants were pulled out from
the tubes and river sand were washed by water. Then, the roots of Z. macrostachya were soaked in 710 mM NaCl solution (treatment) and
distilled water (control), respectively. Three biological replications
were set for each treatment, three treatments were marked as treatment1,
treatment2, treatment3; and three controls were marked as
control1, control2, control3. After 8 h, 3 cm
root tip were harvested and stored in liquid nitrogen for extracting RNA. Total RNA was extracted from three
treatment and three control samples, respectively, using a mirVana miRNA
Isolation Kit (Ambion) in accordance with the manufacturer’s protocol.
Preparation and sequencing of cDNA library
The integrity of RNA was assessed by the Agilent 2100
Bioanalyzer (Agilent Technologies, Santa Clara, CA, USA). Afterwards, all
samples with the RNA Integrity Number (RIN) of ≥7 were used for
constructing cDNA libraries. Six cDNA libraries were constructed by the TruSeq
Stranded mRNA LTSample Prep Kit (Illumina, San Diego, CA, USA) in accordance
with the manufacturer’s protocol. The above libraries were sequenced using the Illumina HiSeq XTen system to generate the
125/150 bp reads with paired-end.
Quality
control and de novo
assembly
Original data (raw reads) were processed by Trimmomatic (Bolger et al. 2014). Reads that contained
ploy-N and were of low quality were eliminated to obtain clean reads. Later,
the adaptor and sequences of low quality were removed, clean reads were
assembled into expression sequence tag clusters (contigs), then de novo assembled into transcript
through Trinity (version: trinityrnaseq_r20131110) according to the paired-end
approach (Grabherr et al. 2011).
Subsequently, transcript with the greatest length was selected as the unigene
based on the sequence length and similarity.
Function
annotations
Unigene functions were annotated through aligning them with SwissProt
protein, NCBI non-redundant protein (NR), Clusters of orthologous groups for
eukaryotic complete genomes (KOG) and Pfam databases
using Blastx (Altschul et al. 1990),
with the cut-off E-value of 10−5. Typically, those proteins
with the greatest hits to the above unigenes were utilized for assigning the
functional annotations. According to SwissProt annotations, the gene ontology
(GO) analysis was performed based on mapping relation between SwissProt and GO
terms. All unigenes were mapped to the Kyoto Encyclopedia of Genes and Genomes
(KEGG) database for annotating the underlying metabolic pathways (Kanehisa et al. 2008).
Unigenes
quantification, analysis of DEGs, cluster analysis, GO and KEGG enrichment analyses
FPKM value (Trapnell et al.
2010), together with the read counts for every unigene, were computed through
eXpress (Roberts and Lior 2013) and bowtie2 (Langmead and Salzberg 2012).
Afterwards, DEGs were identified through nbinom Test and DESeq (Anders and
Huber 2012) functions estimate Size
Factors. Unigenes with P-value of ≤0.05 and log2foldChange of ≥1 were selected as thresholds of significant
DEGs. Later, DEGs were performed hierarchical cluster analysis for exploring
the expression profiles of transcripts. Then, the DEGs were carried out for GO
enrichment and KEGG pathway analyses using the R software based on the
hypergeometric distribution.
Annotation of TFs
The Hidden Markov Model (HMM) motif sequences of TFs were
obtained based on Plant Transcription Factor Database (TFDB), which contained
58 TF families of green plant (Jin et al.
2014). Specifically, one unigene with ≥90% sequence homology was
annotated as the putative TF (E-value of ≤1E-10). Thereafter,
DEGs, together with the candidate TFs in reference transcriptome, were
clustered according to the TF families.
Quantitative
real-time PCR (qRT-PCR)
For further validating the results of transcriptomic
analyses,
10 DEGs including 1 POD-related, 2 SOD-related, 2 KEA-related, 2
KUP-related, 1 KOC-related, 1 Proline-related and 1 Betaine-related DEGs were
selected for
qRT-PCR analysis.
Z. macrostachya were in exposure of 710 mM NaCl solution and distilled water
for 8 h, respectively, then the root RNA were extracted according to the method mentioned above.
Reverse transcription reactions were performed using SuperScript III
reverse transcriptase (Invitrogen, Grand Island, NY, USA) following the
manufacturer's instructions. Primers (Table 1) for qRT-PCR were
designed using the Premier v5.0 software (Premier Biosoft, Palo Alto, CA, USA)
with β-actin genes as the internal controls. The qRT-PCR was carried out
by the Two-color Real-time PCR Detection System (Bio-Rad, USA) following amplification protocol:
3 min at 95°C, and then 3 s at 95°C and 30 s at 60°C for 45 cycles. All
reactions were performed in triplicate, the relative expression levels of the
selected unigenes normalized to β-actin was calculated using the 2
−△△Ct method
(Livak and Schmittgen 2001).
%20591-600_files/image002.jpg)
Fig. 1: Unigene
length distribution
Results
Sequencing
results and assembly
A mean of 49.49 million raw reads were produced from
controls, whereas 49.44 million raw reads were produced from 710 mM NaCl treatments (Table 2). Over
93.99% of all raw reads possessed the Phred-like quality score of Q30 level
(error probability=1‰). An average of 48.23 million (controls) and 48.19
million (treatments) clean reads were obtained, and the valid base ratio and GC
content were above 95.28 and 51.05%, respectively, indicating high quality of
sequencing and cDNA library establishment. A total of 4,7325 unigenes with mean
length of 1397 bp and the N50 of 2169 bp were obtained. As shown in Fig. 1,
24,192 unigenes had the length of 301–1000 bp, 10,317 had the length of
1001–2000 bp, and 11,431 had the length of >2000 bp.
Clean reads of each sample were
mapped to Z. macrostachya unigenes to
confirm the quality of sequencing for those six samples (Table 3). As observed,
the average total mapped reads proportion was 87.61% (range, 86.98–88.44%), and
the multi-position matched reads accounted for approximately 21.65%, and 64.89–67.06% reads were uniquely matched unigenes in all the
six samples. An average of 81.43% reads was mapped in pairs.
Functional
annotations of unigenes
Unigenes were annotated by aligning to the publicly
accessible databases (Table 4). Estimated number of 32,033 (67.69%) unigenes
were aligned in the NR database, 24,155 (51.04%) in the Swissprot database,
11,601 (24.51%) in the KEGG database, 17,906 (26.78%) in the KOG database,
22,309 (47.14%) in the GO database, and 51(0.11%) in the Pfam database.
Altogether 32,542 unigenes (68.76%) were annotated against at least one of the
following databases, including NR, SwissProt, KEGG, Pfam, GO and KOG.
According to the homology of
sequences, 22,309 unigenes were clustered to 3 major GO classifications,
including biological process, cellular component and molecular function (Fig.
2). These unigenes were subdivided into 57 GO terms. Among biological process
classifications, the term of “biological regulation”, “cellular process”,
“metabolic process”, “regulation of biological process” and “responses to
stimulus” were the dominant clusters, whereas only a few unigenes belonged to
the “biological adhesion”, “cell killing”, “locomotion” and “rhythmic process”
term. With regard to the cellular component classification, the term of “cell”,
“cell part” and “organelle” had the bigger unigenes proportions. As to the
molecular function, most unigenes were classified into the “binding”,
“catalytic activity” and “transporter activity” term. Interestingly, 225
unigenes were categorized into the “antioxidant activity” term, while 894 as
the “transporter activity” term.
Table 1: The qRT-PCR primers for detecting the accuracy of transcriptome data
|
Unigene number |
Unigene name |
Annotated as |
Forward primers (5'-3') |
Reverse primers (5'-3') |
|
1 |
TRINITY_DN11437_c0_g1_i1_2 |
POD |
GCAAGTCAAGCACCTCAACCT |
TGATTGGAATGGTCTGCTGGGA |
|
2 |
TRINITY_DN10853_c0_g1_i2_1 |
SOD |
GACTGAATCTCTACGCCTGT |
GACGACGTATGGCACCAGAG |
|
3 |
TRINITY_DN19735_c0_g1_i6_2 |
SOD |
AGCTGCTCCATTGCCATTCCT |
CCGAACCTCCTGTAAGTCAACC |
|
4 |
TRINITY_DN19255_c0_g2_i4_1 |
KEA |
TAGTAAGGGAGAATCTTTGCAAG |
ATGGCAGATCCGTGACAAGTC |
|
5 |
TRINITY_DN22016_c0_g1_i24_2 |
KEA |
TCACCGCCATCCCAGTCATC |
CCCAAATACACGCACTCCTG |
|
6 |
TRINITY_DN16474_c0_g3_i4_2 |
KUP |
CTTGTCGACCAACATTCACACC |
CGAGTAAGCAACGGTCCA |
|
7 |
TRINITY_DN16569_c0_g2_i2_1 |
KUP |
ACTGGACGAGGAGCAACACT |
TCTCACTCACTCCTCCAGAGC |
|
8 |
TRINITY_DN16911_c0_g2_i5_2 |
KOC |
TCCCAAACCAGCTTCACCGATT |
CACCAGAGTATGCCTCGACAT |
|
9 |
TRINITY_DN20583_c0_g1_i5_1 |
Free proline |
TTACAACGGTTCGCAACTGT |
TGCGTGGACTACTCAGAGACC |
|
10 |
TRINITY_DN15030_c0_g1_i1_1 |
Glycine betaine |
TCTCCTCATCCTCACCTCAT |
ATTGCCAGTTCCTAGTGTTCC |
Table 2: Summary of
the sequence analysis
|
Sample |
Raw reads |
Raw bases |
Clean reads |
Clean bases |
Q30 |
GC |
|
|
7430464200 |
48215536 |
7090647459 |
95.43% |
94.42% |
51.16% |
||
|
Control2 |
49533550 |
7430032500 |
48398632 |
7115121942 |
95.76% |
94.68% |
51.84% |
|
Control3 |
49424632 |
7413694800 |
48083930 |
7063618562 |
95.28% |
93.99% |
51.40% |
|
7452987300 |
48368100 |
7120357990 |
95.54% |
94.22% |
51.59% |
||
|
Treatment2 |
49108932 |
7366339800 |
47923128 |
7058090438 |
95.82% |
94.58% |
51.05% |
|
Treatment3 |
49536672 |
7430500800 |
48291630 |
7112106699 |
95.72% |
94.38% |
51.26% |
Table 3: Mapping of clean reads to unigenes
|
Terms\Samples |
Control1 |
Control2 |
Control3 |
Treatment1 |
Treatment2 |
Treatment3 |
|
Total reads |
48215536 (100.00%) |
48398632 (100.00%) |
48083930 (100.00%) |
48368100 (100.00%) |
47923128 (100.00%) |
48291630 (100.00%) |
|
Total mapped reads |
42448603 (88.04%) |
42802853 (88.44%) |
42159835 (87.68%) |
42069655 (86.98%) |
41718189 (87.05%) |
42224960 (87.44%) |
|
Multiple mapped |
10306193 (21.38%) |
10344313 (21.37%) |
10207699 (21.23%) |
10490280 (21.69%) |
10622363 (22.17%) |
10640511 (22.03%) |
|
Uniquely mapped |
32142410 (66.66%) |
32458540 (67.06%) |
31952136 (66.45%) |
31579375 (65.29%) |
31095826 (64.89%) |
31584449 (65.40%) |
|
Reads mapped in proper pairs |
39451806 (81.82%) |
39954228 (82.55%) |
39326280 (81.79%) |
38891334 (80.41%) |
38763888 (80.89%) |
39160458 (81.09%) |
Table 4: Blast
analysis of non-redundant unigenes against public databases
|
Databases for annotation |
Unigene numbers |
Percentages |
300≤ Length<1000nt |
Length≥1000nt |
|
Annotated in NR |
32033 |
67.69 % |
10665 |
21368 |
|
Annotated in Swissprot |
24155 |
51.04 % |
6615 |
17540 |
|
Annotated in KEGG |
11601 |
24.51 % |
3571 |
8030 |
|
Annotated in KOG |
17906 |
37.84 % |
4774 |
13132 |
|
Annotated in GO |
22309 |
47.14 % |
6231 |
16078 |
|
Annotated in Pfam |
51 |
0.11 % |
41 |
10 |
|
Annotated in at least one
database |
32542 |
68.76% |
21547 |
12995 |
%20591-600_files/image004.png)
Fig. 3: KOG
classification of unigenes
The KOG database was searched to
identify unigenes to predict and classify their functions. Based on the
sequence homology, 17,906 sequences showed a KOG classification (Fig. 3). Among
the 25 KOG clusters, the cluster of “general function prediction only” (3674,
20.52%) was the largest group, followed by the “posttranslational modification,
protein turnover, chaperones” (2180, 12.17%), the “signal transduction
mechanisms” (1821, 10.17%) and the “translation, ribosomal structure and
biogenesis” (1266, 7.07%). By contrast, the clusters of the “cell motility” (3,
0.17‰), the “nuclear structure” (51, 0.28%) and the “extracellular structures”
(52, 0.29%) were smaller groups.
A total of 11,601 annotated unigenes were mapped to the reference
canonical pathways in the KEGG database (Fig. 4). The most represented pathways
were “metabolism” (5711 unigenes), “genetic information processing” (2822
unigenes), “cellular processes” (1430 unigenes) and “environmental information
processing” (1298 unigenes). In addition, “translation” (1335 unigenes),
“signal transduction” (1216 unigenes) and “carbohydrate metabolism” (1049
unigenes) were the pathways most closely related to the mapped unigenes, which
suggested the presence of active pathways in Z. macrostachya.
%20591-600_files/image006.png)
Fig. 4: KEGG
categorization of unigenes
%20591-600_files/image008.jpg)
Fig. 5: Volcano
plots of DEGs
Recognition and annotation of DEGs
A total of 14,558 DEGs were recognized, among which, 7972 were
up-regulated and 6586 were down-regulated (Fig. 5). According to our
findings, 7424 DEGs were enriched into 64 GO terms. With regard to biological
process, the term of “cellular process”, “metabolic process”, “biological
regulation”, “regulation of biological %20591-600_files/image010.jpg)
Fig. 6:
GO categorization of DEGs in Z. macrostachya. under salt stress
%20591-600_files/image012.jpg)
Fig. 2: GO
categorization of unigenes
process” and “response to stimulus” had the greater
unigenes proportions. As to cellular component, the “membrane”, the “cell part”
and the “cell” term were the dominant clusters. The “binding”, “catalytic
activity” and “transporter activity” term contained more unigenes in the
molecular function category (Fig. 6).
In
addition, DEGs were performed for KEGG analysis. The results suggested that,
2093 DEGs were assigned %20591-600_files/image014.png)
Fig. 7: KEGG
categorization of DEGs in Z. macrostachya
under salt stress
with the KEGG ID, which were
classified into 197 pathways (Fig. 7). Typically, the most representative
pathways were the “environmental information processing-signal transduction”
(370 unigenes), the “metabolism-amino acid metabolism” (248 unigenes), the
“genetic information processing-translation” (224 unigenes), the
“metabolism-global and overview maps” (223 unigenes) and the “metabolism-lipid
metabolism” (219 unigenes).
DEGs
associated with ROS scavenging
A total of 24 DEGs were identified to encode ROS
scavenging-associated enzymes (Table 5). Among which, 12 encoded PODs and
constituted the largest group. The expression of TRINITY_DN11379_c0_g1_i1_2,
TRINITY_DN11437_c0_g1_i1_2 and TRINITY_DN12963_c0_g1_i2_2 increased under salt
stress, with high transcript levels. Four genes encoded SODs, the expression of
one gene increased, while that of the other three genes decreased. Two genes
encoded CATs, among which, TRINITY_DN21281_c0_g4_i5_2 showed higher transcript
abundance, but TRINITY_DN21281_c0_g1_i8_2 was down-regulated under salt
stress. Additionally, five genes were found to encode Ascorbate peroxidase
(APX) and one gene encoded glutathione peroxidase (GPX). The encoded APX gene
TRINITY_DN19584_c0_g2_i1_2 was up-regulated by nearly 26 times (FC 26.80), with
the highest transcript levels.
DEGs
involved in ion transportation
Fifty-two DEGs were identified to be related to the
regulation of K+ transportation (Table 6), which accounted for the
greatest proportion of identified genes. 37 out of them were up-regulated, while
15 were down-regulated under salt stress. Eleven of them were involved in
K+ efflux antiporter (KEA), twenty were related to K+
transmembrane transporter (KUP), three were associated with K+
channel (AKT), nine were involved in the outward rectifying K+
channel (KOC) and nine were related to cyclic nucleotide-gated channel (CNGC).
Three genes were identified to encode the plasma membrane P-ATPases (PM-H+-ATPases),
while five encoded the vacuolar V-ATPases (V-H+-ATPases), and these
unigenes were up-regulated. Two unigenes that encoded Na+/H+
exchanger (NHX) were also identified, but their expressions levels were
down-regulated.
DEGs
associated with osmotic adjustment
In this study, 12 DEGs were identified to involve in free
proline metabolism, eight out of them were up-regulated, while four were
down-regulated under salt stress (Table 7). Four DEGs were identified to
involve in Glycine betaine metabolism, 3 out of them were up-regulated, only
TRINITY_DN15030_c0_g1_i1_1 was down-regulated under salt stress (Table 7).
DEGs
related to TFs
Table 5: DEGs related
to the ROS scavenging system
|
Gene ID |
Log2FC |
P-value |
Gene ID |
Log2FC |
P-value |
|
POD |
|
|
SOD |
|
|
|
TRINITY_DN10724_c0_g1_i2_1 |
-1.61 |
1.46E-02 |
TRINITY_DN10853_c0_g1_i2_1 |
-2.67 |
1.13E-11 |
|
TRINITY_DN11379_c0_g1_i1_2 |
4.18 |
2.69E-11 |
TRINITY_DN19735_c0_g1_i6_2 |
1.97 |
3.98E-04 |
|
TRINITY_DN11437_c0_g1_i1_2 |
5.43 |
9.90E-03 |
TRINITY_DN699_c0_g1_i1_2 |
-1.13 |
5.85E-06 |
|
TRINITY_DN11915_c0_g1_i1_1 |
-2.11 |
1.28E-25 |
TRINITY_DN7158_c0_g1_i2_1 |
-2.22 |
1.73E-11 |
|
TRINITY_DN12563_c0_g1_i1_1 |
-3.11 |
2.07E-16 |
APX |
|
|
|
TRINITY_DN12830_c0_g1_i3_1 |
1.20 |
4.15E-07 |
TRINITY_DN18445_c0_g1_i11_1 |
-4.10 |
1.56E-08 |
|
TRINITY_DN12963_c0_g1_i2_2 |
4.89 |
5.57E-14 |
TRINITY_DN18597_c0_g2_i1_2 |
1.56 |
1.30E-02 |
|
TRINITY_DN13326_c0_g1_i1_2 |
1.42 |
2.82E-02 |
TRINITY_DN19584_c0_g1_i3_2 |
1.99 |
1.39E-29 |
|
TRINITY_DN13354_c0_g1_i1_1 |
-3.66 |
5.97E-12 |
TRINITY_DN19584_c0_g2_i1_2 |
4.74 |
1.61E-03 |
|
TRINITY_DN14879_c0_g3_i1_1 |
-2.59 |
1.10E-15 |
TRINITY_DN9766_c0_g1_i2_2 |
-2.91 |
7.99E-46 |
|
TRINITY_DN14924_c0_g1_i1_1 |
-2.67 |
2.84E-02 |
GPX |
|
|
|
TRINITY_DN15006_c1_g1_i1_1 |
2.78 |
7.16E-11 |
TRINITY_DN13470_c0_g1_i1_1 |
1.21 |
9.32E-05 |
|
CAT |
|
|
|
|
|
|
TRINITY_DN21281_c0_g1_i8_2 |
-2.45 |
5.95E-44 |
|
|
|
|
TRINITY_DN21281_c0_g4_i5_2 |
2.28 |
3.62E-32 |
|
|
|
Table 6: DEGs related to K+ and Na+ ion transport
|
Gene ID |
Log2FC |
P-value |
Gene ID |
Log2FC |
P-value |
|
KEA |
|
|
KOC |
|
|
|
TRINITY_DN14177_c0_g1_i4_2 |
-1.65 |
9.15E-22 |
TRINITY_DN14314_c0_g1_i7_1 |
1.52 |
1.39E-21 |
|
TRINITY_DN14664_c0_g1_i3_2 |
1.82 |
2.80E-07 |
TRINITY_DN16211_c0_g1_i2_2 |
2.52 |
5.34E-06 |
|
TRINITY_DN15615_c0_g1_i2_2 |
-1.10 |
2.96E-04 |
TRINITY_DN16911_c0_g2_i5_2 |
7.38 |
1.02E-203 |
|
TRINITY_DN17765_c1_g1_i9_1 |
3.65 |
5.91E-08 |
TRINITY_DN19502_c0_g3_i1_2 |
1.21 |
5.76E-06 |
|
TRINITY_DN17765_c1_g2_i5_1 |
-1.06 |
3.66E-02 |
TRINITY_DN20260_c0_g1_i20_2 |
2.08 |
3.86E-32 |
|
TRINITY_DN19255_c0_g1_i3_1 |
-3.11 |
3.36E-06 |
TRINITY_DN20505_c0_g1_i1_1 |
3.24 |
4.61E-54 |
|
TRINITY_DN19255_c0_g2_i4_1 |
-3.14 |
7.72E-05 |
TRINITY_DN20585_c1_g2_i17_2 |
1.17 |
2.61E-05 |
|
TRINITY_DN21559_c1_g1_i5_2 |
2.23 |
1.58E-18 |
TRINITY_DN22190_c2_g1_i9_2 |
3.48 |
6.03E-49 |
|
TRINITY_DN21559_c1_g5_i3_2 |
1.85 |
1.67E-05 |
TRINITY_DN22190_c2_g2_i1_2 |
3.45 |
3.48E-15 |
|
TRINITY_DN22016_c0_g1_i24_2 |
3.69 |
2.89E-12 |
CNGC |
|
|
|
TRINITY_DN22016_c0_g2_i3_2 |
1.12 |
1.64E-04 |
TRINITY_DN12577_c0_g1_i1_2 |
2.10 |
6.61E-23 |
|
KUP |
|
|
TRINITY_DN19156_c0_g1_i14_2 |
1.99 |
1.11E-05 |
|
TRINITY_DN15916_c0_g1_i8_1 |
1.80 |
5.39E-04 |
TRINITY_DN19236_c0_g1_i7_2 |
2.09 |
1.67E-37 |
|
TRINITY_DN16437_c0_g1_i1_2 |
1.74 |
1.27E-07 |
TRINITY_DN20366_c0_g3_i7_1 |
-1.15 |
1.52E-04 |
|
TRINITY_DN16474_c0_g2_i1_2 |
2.08 |
6.62E-05 |
TRINITY_DN20693_c0_g1_i17_1 |
-2.42 |
3.16E-08 |
|
TRINITY_DN16474_c0_g3_i4_2 |
4.16 |
9.84E-103 |
TRINITY_DN21311_c0_g1_i4_2 |
1.53 |
7.62E-09 |
|
TRINITY_DN16569_c0_g2_i2_1 |
-2.69 |
1.11E-08 |
TRINITY_DN21317_c0_g1_i1_2 |
1.34 |
1.39E-09 |
|
TRINITY_DN16999_c0_g1_i4_1 |
2.05 |
8.74E-26 |
TRINITY_DN21317_c0_g2_i2_2 |
1.27 |
1.29E-07 |
|
TRINITY_DN18621_c0_g1_i6_2 |
1.78 |
4.17E-09 |
TRINITY_DN22290_c0_g1_i1_2 |
1.58 |
1.26E-06 |
|
TRINITY_DN18799_c0_g1_i2_2 |
-1.65 |
6.76E-07 |
P-ATPase |
|
|
|
TRINITY_DN19410_c0_g1_i7_1 |
1.17 |
1.81E-02 |
TRINITY_DN20505_c0_g2_i2_1 |
4.50 |
1.70E-135 |
|
TRINITY_DN19574_c0_g1_i1_2 |
2.07 |
3.17E-14 |
TRINITY_DN20585_c1_g2_i17_2 |
1.17 |
2.61E-05 |
|
TRINITY_DN19850_c0_g1_i9_2 |
2.38 |
3.45E-23 |
TRINITY_DN22190_c2_g1_i9_2 |
3.48 |
6.03E-49 |
|
TRINITY_DN19858_c0_g1_i37_2 |
2.97 |
4.48E-65 |
V-ATPase |
|
|
|
TRINITY_DN19995_c2_g1_i15_1 |
-1.30 |
2.63E-06 |
TRINITY_DN13353_c0_g1_i2_2 |
5.91 |
2.81E-35 |
|
TRINITY_DN21041_c0_g1_i3_2 |
1.94 |
5.25E-31 |
TRINITY_DN15633_c0_g1_i2_2 |
1.64 |
2.75E-09 |
|
TRINITY_DN22105_c0_g1_i4_2 |
3.27 |
1.64E-54 |
TRINITY_DN14314_c0_g1_i7_1 |
1.52 |
1.39E-21 |
|
TRINITY_DN24434_c0_g1_i1_1 |
-2.95 |
1.47E-04 |
TRINITY_DN20260_c0_g1_i20_2 |
2.08 |
3.86E-32 |
|
TRINITY_DN4618_c0_g1_i1_1 |
-3.05 |
6.68E-35 |
TRINITY_DN20505_c0_g1_i1_1 |
3.24 |
4.61E-54 |
|
TRINITY_DN5960_c0_g1_i1_1 |
1.95 |
4.30E-20 |
NHX |
|
|
|
TRINITY_DN7740_c0_g1_i1_1 |
-3.82 |
6.23E-05 |
TRINITY_DN14878_c0_g1_i1_2 |
-1.02 |
1.16E-03 |
|
TRINITY_DN7928_c0_g1_i2_1 |
-2.28 |
1.36E-10 |
TRINITY_DN9153_c0_g1_i2_1 |
-2.51 |
2.12E-04 |
|
AKT |
|
|
|
|
|
|
TRINITY_DN13806_c0_g1_i1_1 |
-2.04 |
3.70E-08 |
|
|
|
|
TRINITY_DN13806_c0_g2_i4_1 |
-2.04 |
1.73E-02 |
|
|
|
|
TRINITY_DN24232_c0_g1_i1_1 |
1.71 |
2.59E-09 |
|
|
|
A total of 2,327 DEGs were annotated to the TFs database,
which belonged to 57 families (Fig. 8). 604 out of them
were annotated to the BHLH family, which accounted for the largest group,
followed by NAC (470 DEGs), MYB-related (402 DEGs) and ERF (356 DEGs) family. The
LSD (5 DFGs), LFY (5 DEGs) and SAP (1 DEGs) family had the least corresponding
unigenes.
Validation of DEGs relative expression by qRT-PCR
Table 7: DEGs related to osmotic
adjustment substances
|
Gene ID |
Log2FC |
P-value |
|
Free proline |
|
|
|
TRINITY_DN12426_c0_g2_i1_1 |
-2.54 |
2.81E-02 |
|
TRINITY_DN16249_c0_g1_i5_2 |
5.09 |
7.53E-51 |
|
TRINITY_DN16737_c0_g1_i1_2 |
4.69 |
3.33E-17 |
|
TRINITY_DN20583_c0_g1_i5_1 |
6.72 |
7.84E-272 |
|
TRINITY_DN22133_c0_g1_i14_2 |
3.78 |
3.74E-60 |
|
TRINITY_DN22200_c0_g1_i19_2 |
7.87 |
0.00E+00 |
|
TRINITY_DN20075_c0_g1_i5_1 |
2.19 |
8.55E-36 |
|
TRINITY_DN20857_c1_g4_i6_1 |
2.40 |
1.83E-47 |
|
TRINITY_DN19370_c0_g1_i3_1 |
-2.85 |
1.21E-08 |
|
TRINITY_DN20069_c0_g1_i5_2 |
1.74 |
8.26E-22 |
|
TRINITY_DN15955_c0_g1_i5_1 |
-1.31 |
5.68E-03 |
|
TRINITY_DN19808_c0_g2_i7_1 |
-2.45 |
7.36E-10 |
|
Glycine betaine |
|
|
|
TRINITY_DN12892_c0_g1_i2_1 |
2.29 |
1.74E-22 |
|
TRINITY_DN15030_c0_g1_i1_1 |
-2.86 |
5.60E-05 |
|
TRINITY_DN18859_c0_g1_i5_2 |
1.82 |
9.97E-04 |
|
TRINITY_DN12184_c0_g1_i3_2 |
2.14 |
9.46E-19 |
%20591-600_files/image016.png)
Fig. 8: Distribution
of TFs family
%20591-600_files/image018.png)
Fig. 9: Expression
pattern validation of selected genes by qRT-PCR
The qRT-PCR
results showed the expression patterns of 10 selected DEGs in good agreement
with the relative expression results at transcriptome level (Fig. 9). This
indicated that transcriptome data of Z.
macrostachya roots obtained
in this experiment were accurate.
Discussion
In this
study, the roots informative transcriptome dataset for Z. macrostachya after NaCl treatment revealed that 68.76% of the
4,7325 unigenes were annotated by BLAST analysis, which suggested that the
sequences of the Z. macrostachya
unigenes generated in the present study were assembled and annotated correctly. A total of 14,558 DEGs were identified from Z. macrostachya roots under salt stress,
these unigenes provided a comprehensive understanding towards the genes
transcription profiles of Z. macrostachya,
and laid a solid foundation for further study of salt-tolerance mechanisms and
identification of new genes in this species.
ROS accumulation induce cytoplasmic
membrane damage (Xiong et al. 2020). Study has
shown that salt-tolerance of plant was related to scavenging capacity of
antioxidative proteins to some extent (Zhang et al. 2012; Bose et al.
2014). SOD is the first-line of defense in resistance to
oxidative injury, which can dismutate •O2- into O2
and H2O2 (Qu et al. 2010), CAT and POD can dismutate H2O2
into H2O and O2 (Xu et
al. 2013), whereas, APX exerts a vital part in the catalysis of H2O2
conversion to H2O, and H2O2 can also be
reduced by GPX in the GPX pathway (Gu et
al. 2018). In this study, these DEGs encodesd enzymatic antioxidants might
be the vital factors involved in ROS elimination in Z. macrostachya.
Ions transport is a crucial factor in plant response to
salt stress. To deal with salt stress, plants have evolved
specific mechanisms for coordinating the Na+ discharge and K+
uptake processes (Guo et al.
2016; Lv et al. 2018). KUP, KEA and
AKT constitute the K+ uptake system in plant, PM-H+-ATPase
and V-H+-ATPase can generate the driving force to extrude Na+
out of cell and compartmentalize Na+ to vacuole, and these are
regarded as the mechanisms that a plant employs to restore the homeostasis
of ions in cell (Chinnusamy et
al. 2006). In our study, 37 up-regulated DEGs involved in K+
uptake, 7 up-regulated DEGs involved in Na+ extrusion was
identified, respectively, indicating these DEGs played critical roles in
reestablishing the cellular K+ and Na+ homeostasis. In
addition, NHX were probably responsible for Na+ sequestration, but 2
DEGs encoding NHX were down-regulated, revealing that
NHXs were not likely involved in Na+ discharge.
It is well-known that plants
synthesize some osmolytes, such as free proline and Glycine betaine, so as to
alleviate the negative effects of toxic ions and osmotic imbalance (Ingram and
Bartels 1996; Ashraf and Foolad 2007). In this study, 11 DEGs involved in free
proline and Glycine betaine metabolism were identified, their expression were
up-regulated, it suggested that these DEGs exerted vital parts in regulating
osmotic balance in Z. macrostachya.
The TFs can modulate the downstream
genes expression, which play important roles in plant response to abiotic stresses
(Capella et al. 2015). Specifically,
some TFs families have been identified from halophytes through RNA-Seq, such as
Iris lactea var. chinensis (Gu et al. 2018), Reaumuria trigyna (Wang et al.
2014) and Suaeda fruticosa (Diray-Arce
et al. 2015). In our study, some DEGs
were annotated to the BHLH, NAC, MYB-relate or ERF TFs
family. It is reported that the expression of TabHLH13 isolated from wheat was
rapidly up-regulated after salt stress (Kim and Kim 2006), while that of bHLH92
in Arabidopsis thaliana was triggered
under drought and high salinity (Jiang and Deyholos 2006; Liu et al. 2014), indicating that the BHLH
TFs family exert vital parts in the salt-tolerance of plant. Previous studies
showed that, the over-expression of NAC TF GhATAF1 improved the salt-tolerance
of cotton (He et al. 2016); besides,
the ERFs TFs family can modulate the responses to abiotic stress and ethylene,
as well as disease resistance, which were achieved through targeting the
downstream promoter GCC-box (Ohme-Takagi and Shinshi 1995; Stockinger et al. 1997). By contrast, the
expression of MYB-related TF AtMYBL in A.
thaliana is activated by salt stress (Zhang et al. 2011). Thus, the TFs identified in our study laid the basis
for more intensive studies.
Conclusion
A total of
4,7325 unigenes were assembled, and 68.76% of the 4,7325 unigenes were
annotated by BLAST analysis. Altogether 14,558 DEGs were identified between
salt stress and control samples, many DEGs were identified to encode enzymes
related to the ROS scavenging system, K+ and Na+
transport proteins, osmotic adjustment substances and TFs, which potentially
played vital roles in salt stress response in Z. macrostachya. Our findings provide useful insights for
salt-tolerance molecular mechanisms of Z.
macrostachya and provide basis for future cloning of salt-tolerance genes
in this species.
Acknowledgements
This work
was funded by the Jiangsu
Key Laboratory for Bioresources of Saline Soils Fund (JKLBS2017009), and
the Six Talent Peaks Project in Jiangsu province (NY-096) in 2017.
Author Contributions
HG Hu and
ZM Zhang designed the experiments; HG Hu, ZM Zhang and SN Min performed the
experiments; HG Hu and ZM Zhang analyzed RNA-Seq data; HG Hu, ZM Zhang and SN
Min wrote the manuscript.
References
Altschul SF, W Gish, W Miller, EW Myers, DJ Lipman
(1990). Basic local alignment search tool.
J Mol Biol 215:403‒410
Anders S, W Huber
(2012). Differential expression of RNA-Seq data at the gene level–the DESeq
package. EMBL 3:1‒23
Ashraf M, MR Foolad (2007). Roles of glycine bataine and
proline in improving plant abiotic stress resistance. Environ Exp Bot 59:206‒216
Bolger AM, M Lohse, B Usadel (2014). Trimmomatic: A flexible trimmer for Illumina sequence data. Bioinformatics 30:2114‒2120
Bose J, A Rodrigo-Moreno,
S Shabala (2014). ROS homeostasis in halophytes in the context of salinity
stress tolerance. J Exp Bot 65:1241‒1257
Capella M, PA Ribone, AL Arce, RL Chan (2015). Arabidopsis thaliana HomeoBox 1(AtHB1), a Homedomain-Leucine Zipper I
(HD-Zip I) transcription factor, is regulated by PHYTOCHROME-INTERACTING FACTOR
1 to promote hypocotyl elongation. New Phytol 207:669‒682
Chinnusamy V, JH Zhu, JK Zhu (2006). Salt stress
signaling and mechanisms of plant salt tolerance. Genet Eng 27:141‒177
Deinlein U, AB Stephan, T Horie, W Luo, G Xu, JI
Schroeder (2014). Plant salt-tolerance mechanisms. Trends Plant Sci 19:371‒379
Diray-Arce J, M Clement, B Gul, MA Khan, BL Nielsen
(2015). Transcriptome assembly, profiling and differential gene expression
analysis of the halophyte Suaeda
fruticosa provides insights into salt tolerance. BMC Genomics 16; Article 353
Flowers TJ, TD Colmer (2008). Salinity tolerance in
halophytes. New Phytol 179:945‒963
Grabherr MG, BJ Haas, M Yassour, JZ Levin (2011). Trinity: Reconstructing a full-length
transcriptome without a genome from RNA-Seq data. Nat Biotechnol 29:644‒652
Gu C, S Xu, Z Wang, L Liu, Y Zhang, Y Deng, S Huang
(2018). De novo sequencing, assembly,
and analysis of Iris lactea var. chinensis roots’ transcriptome in
response to salt stress. Plant Physiol
Biochem 125:1‒12
Guo HL, Y Wang, DD Li, JB Chen, JQ Zong, ZY Wang, X Chen, JX Liu (2016). Growth response and ion regulation of seashore paspalum
accessions to increasing salinity. Environ Exp Bot 131:137‒145
He X, L Zhu,
L Xu, W Guo, X Zhang (2016). GhATAF1, a NAC transcription factor, confers
abiotic and biotic stress responses by regulating phytohormonal signaling
Networks. Plant Cell Rep 35:2167‒2179
Hu HG, ZM Zhang (2009). Study on the physiologicai
adaptability of Zoysia macrostachya to
salt stress. J Anhui Agric Sci
37:6413‒6415
Hu HG, ZM Zhang (2010). Study on responses to salt stress
and critical salt concentration of Zoysia macrostachya. Northern Hortic 3:80‒83
Hu HG, ZM Zhang, TX Sun, YY Han (2016). Effects of mixed
saline stress on the absorption and transportation of K+ and Na+
of Zoysia Willd. Chin J Trop Crops 37:2145‒2149
Ingram J, D Bartels (1996). The molecular basis of
dehydration tolerance in plants. Annu
Rev Plant Biol 7:377‒403
Jiang Y, MK Deyholos (2006).
Comprehensive transcriptional profiling of NaCl-stressed Arabidopsis roots reveals novel classes of responsive genes. BMC Plant Biol 6:25‒30
Jin J, H Zhang, L Kong, G Gao, J Luo (2014). PlantTFDB
3.0: A portal for the functional and evolutionary study of plant transcription
factors. Nucl Acids Res 42:1182‒1187
Kanehisa M, M Araki, S Goto, M
Hattori (2008). KEGG for linking genomes to life and the environment. Nucl Acids Res 36:480‒484
Kim J, HY Kim (2006). Functional analysis of a
calcium-binding transcription factor involved in plant salt stress signaling. FEBS Lett 580:5251‒5256
Langmead B, SL Salzberg (2012). Fast gapped-read alignment with Bowtie 2.
Nat Meth 9:357‒359
Liu W, H Tai, S Li, W Gao, M Zhao, CX Xie, W Li (2014).
bHLH122 is important for drought and osmotic stress resistance in Arabidopsis and in the repression of ABA
catabolism. New Phytol 201:1192‒1204
Liu ZX, YJ Hua, SN Wang, XH Liu, LS Zou, CH Chen, H Zhao,
Y Yan (2020). Analysis of the Prunellae spica transcriptome under salt stress. Plant Physiol Biochem 156:314‒322
Livak KJ,
TD Schmittgen (2001). Analysis of relative gene expression data using real-time
quantitative PCR and the 2-ΔΔCT method. Methods 25:402‒408
Lv XY, Y Jin, YG Wang (2018). De novo transcriptome assembly and identification of
salt-responsive genes in sugar beet M14. Comput
Biol Chem 75:1‒10
Meera SP, A Augustine (2020). De novo transcriptome analysis of Rhizophora mucronata Lam.
Furnishes evidence for the existence of glyoxalase system correlated to
glutathione metabolic enzymes and glutathione regulated transporter in salt
tolerant mangroves. Plant Physiol Biochem
155:683‒696
Miller G, V Shulaev, R
Mittler (2008). Reactive oxygen signaling and abiotic stress. Physiol Plantarum 133:481‒489
Ohme-Takagi M, H
Shinshi (1995). Ethylene-inducible DNA binding proteins that interact with an
ethylene-responsive element. Plant Cell
7:173‒182
Roberts A, P Lior (2013). Streaming
fragment assignment for real-time analysis of sequencing experiments. Nat Meth 10:71‒73
Stockinger EJ, SJ Gilmour, MF Thomashow (1997). Arabidopsis thaliana CBF1 encodes an AP2
domain-containing transcriptional activator that binds to the C-repeat/DRE, a
cis-acting DNA regulatory element that stimulates transcription in response to
low temperature and water deficit. Proc
Natl Acad Sci USA 94:1035‒1040
Tang XL, XM
Mu, HB Shao, HY Wang, M Brestic (2015). Global plant-responding mechanisms to salt
stress: Physiological and molecular levels and implications in biotechnology. Crit Rev Biotechnol 35:425‒437
Trapnell C, BA Williams, G Pertea, A Mortazavi, G
Kwan, MJ van Baren, SL Salzberg, BJ Wold, L Pachter (2010). Transcript assembly and quantification by RNA-Seq reveals unannotated
transcripts and isoform switching during cell differentiation. Nat Biotechnol 28:511‒515
Wang J, LL Zheng,
TP Gu, XF Wang, YC Wang (2014). Cloning and expression analysis of two WRKY transcription factors from the
rare recretohalophyte Reaumuria trigyna. Acta Pratacult
Sin 4:122‒129
Xiong X, YQ Wei, JH Chen, N Liu, YJ Zhang (2020). Transcriptome analysis of genes and pathways associated
with salt tolerance in alfalfa under non-uniform salt stress. Plant Physiol
Biochem 151:323‒333
Xu R, M Yamada, H Fujiyama
(2013). Lipid peroxidation and antioxidative enzymes of two turfgrass species
under salinity stress. Pedosphere
23:213‒222
Zhang F, JL Ding, T Tashpolat, BC Wang (2007). Analysis on
characteristics of soil salinization in the arid regions: A case study in the
delta oasis of Weigan and Kuqa rivers. Acta
Pratacult Sin 16:34‒40
Zhang H, B Han, T Wang, SX Chen, HY Li, YH Zhang, SJ Dai
(2012). Mechanisms of plant salt response: Insights from proteomics. J Proteome Res 11:49‒67
Zhang X, HW Ju, MS Chung, P Huang, SJ Ahn, CS Kim (2011).
The R-R-type MYB-like transcription factor, AtMYBL, is involved in promoting
leaf senescence and modulates an abiotic stress response in Arabidopsis. Plant Cell Physiol 52:138‒148
Zhu JK (2001). Plant salt tolerance. Trends Plant Sci 6:66‒71
Zhu JK
(2003). Regulation of ion homeostasis under salt stress. Curr Opin Plant Biol 6:441‒445